NASA/TM-1999-104606, Vol. 15 



Technical Report Series on 

Global Modeling and Data Assimilation 

Max J. Suarez, Editor 


Volume 15 

A Solar Radiation Parameterization for 
Atmospheric Studies 

Ming-Dah Chou and Max J. Suarez 



Climate and Radiation Branch 
NASA Seasonal to Interannual Prediction Project 


June 1999 




The NASA STI Program Office ... in Profile 


Since its founding, NASA has been dedicated to 
the advancement of aeronautics and space 
science. The NASA Scientific and Technical 
Information (STI) Program Office plays a key 
part in helping NASA maintain this important 
role. 

The NASA STI Program Office is operated by 
Langley Research Center, the lead center for 
NASA’s scientific and technical information. 

The NASA STI Program Office provides access 
to the NASA STI Database, the largest collection 
of aeronautical and space science STI in the 
world. The Program Office is also NASA’s 
institutional mechanism for disseminating the 
results of its research and development activi- 
ties. These results are published by NASA in the 
NASA STI Report Series, which includes the 
following report types: 

• TECHNICAL PUBLICATION. Reports of 
completed research or a major significant 
phase of research that present the results of 
NASA programs and include extensive data or 
theoretical analysis. Includes compilations of 
significant scientific and technical data and 
information deemed to be of continuing 
reference value. NASA’s counterpart of 
peer-reviewed formal professional papers but 
has less stringent limitations on manuscript 
length and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 
and technical findings that are preliminary or 
of specialized interest, e.g., quick release 
reports, working papers, and bibliographies 
that contain minimal annotation. Does not 
contain extensive analysis. 

• CONTRACTOR REPORT. Scientific and 
technical findings by NASA-sponsored 
contractors and grantees. 


• CONFERENCE PUBLICATION. Collected 
papers from scientific and technical 
conferences, symposia, seminars, or other 
meetings sponsored or cosponsored by NASA. 

• SPECIAL PUBLICATION. Scientific, techni- 
cal, or historical information from NASA 
programs, projects, and mission, often con- 
cerned with subjects having substantial public 
interest. 

• TECHNICAL TRANSLATION. 
English-language translations of foreign scien- 
tific and technical material pertinent to NASA’s 
mission. 

Specialized services that complement the STI 
Program Office’s diverse offerings include creat- 
ing custom thesauri, building customized data- 
bases, organizing and publishing research results . . . 
even providing videos. 

For more information about the NASA STI Program 
Office, see the following: 

• Access the NASA STI Program Home Page at 
http://www.sti.nasa.gov/STI-homepage.html 

• E-mail your question via the Internet to 
help@sti.nasa.gov 

• Fax your question to the NASA Access Help 
Desk at (301) 621-0134 

• Telephone the NASA Access Help Desk at 
(301) 621-0390 

• Write to: 

NASA Access Help Desk 
NASA Center for AeroSpace Information 
7121 Standard Drive 
Hanover, MD 21076-1320 



NASA/TM-1999-104606, Vol. 15 



Technical Report Series on 

Global Modeling and Data Assimilation 

MaxJ. Suarez, Editor 

Goddard Space Flight Center, Greenbelt, Maryland 


Volume 15 

A Solar Radiation Parameterization for 
Atmospheric Studies 

Ming-Dah Chou 
MaxJ. Suarez 

Goddard Space Flight Center 
Greenbelt, Maryland 


National Aeronautics and 
Space Administration 

Goddard Space Flight Center 

Greenbelt, Maryland 2077 1 


June 1999 



Available from: 


NASA Center for AeroSpace Information 
7121 Standard Drive 
Hanover, MD 21076-1320 
Price Code: A17 

National Technical Information Service 
5285 Port Royal Road 
Springfield, VA 22161 
Price Code: A10 



Abstract 


This NASA Technical Memorandum documents the solar radiative transfer model 
(CLIRAD-SW) developed at the Goddard Climate and Radiation Branch, which has been 
applied to various atmospheric models used in the Goddarad Laboratory for Atmospheres, 
including a general circulation model, a mesoscale model, and a cloud ensemble model. It 
includes the absorption due to water vapor, 0 3 , 0 2 , C0 2 , clouds, and aerosols. Interactions 
among the absorption and scattering by clouds, aerosols, molecules (Rayleigh scattering), and 
the surface are fully taken into account. Fluxes are integrated virtually over the entire 
spectrum, from 0.175 pm to 10 pm. 

Depending upon the nature of absorption, different approaches are applied to different 
absorbers. In the ultraviolet (UV) and photosynthetically active (PAR) region, the spectrum is 
divided into 8 bands, and single 0 3 absorption coefficient and Rayleigh scattering coefficient 
are used for each band. In the infrared, the spectrum is divided into 3 bands, and the k- 
distribution method is applied. Ten absorption coefficients are used in each band. The flux 
reduction due to 0 2 is derived from a simple function, while the flux reduction due to C0 2 is 
derived from precomputed tables. 

Reflection and transmission of a cloud and aerosol-laden layer are computed using the 
8-Eddington approximation. Fluxes are then computed using the two-stream adding 
approximation. For a cloud layer, the optical thickness is parameterized as a function of 
cloud water/ice amount and the effective particle size, whereas the single-scattering albedo and 
asymmetry factor are parameterized as a function of the effective particle size. The 
parameterizations are applied separately to water and ice particles. A maximum-random 
approximation is adopted for the overlapping of clouds at different heights. Aerosol optical 
properties are specified input parameters to the radiation model. 

A special feature of this model is that absorption due to a number of minor absorption 
bands is included. Individually the absorption in those minor bands is small, but collectively 
the effect is large, ~10 % of the atmospheric heating. Integrated over all spectral bands and 
all absorbers, the surface heating is computed accurately to within a few watts per meter 
squared of high spectral-resolution calculations, and the atmospheric heating rate between 
0.01 hPa and the surface is accurate to within 5%. 

A complete computer code of the CLIRAD-SW model with sample calculations is 
available at http://climate.gsfc.nasa.gov/~chou/clirad_sw 
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1. INTRODUCTION 


Solar radiation is the ultimate source of energy for the Earth’s climate. A small error 
in heating rate calculations might induce large errors in model-simulated climate. It is very 
important that solar heating is calculated accurately in an atmospheric model, not only for the 
total heating but also for the horizontal and vertical distributions of the heating. One of the 
major difficulties in the calculation of atmospheric radiative heating is that it requires a huge 
amount of computing time using detailed spectral and angular integrations. There are tens of 
thousands of molecular lines across the spectrum. The half-width of a line is ~0.05 cm 1 in 
the troposphere, which is much narrower than the mean spacing of lines. As a result, the 
absorption coefficient varies rapidly with wavenumbers. Detailed spectral integration of 
radiative transfer requires >10 6 spectral points. When scattering due to clouds and aerosols is 
involved, angular integration becomes the major computational burden. Adding to the 
complexity of the radiative transfer calculations are the overlapping of absorption due to 
multiple absorbers, the overlapping of clouds at various heights, and the wide ranges of cloud 
microphysical and optical properties. 

There are several solar radiation parameterizations with various complexities available 
for use in atmospheric models (Fouquart et al., 1991 and the references therein). At the 
NASA/Goddard Space Flight Center, Climate and Radiation Branch, we have developed a 
model for solar radiative transfer, which can be traced back at least to the past 10 years. 
Depending upon the spectral properties of absorbers and scatterers, as well as the heights 
where heating is important, various parameterizations are applied. This model has been 
applied to the various atmospheric model used in the Goddarad Faboratory for Atmospheres, 
which include a general circulation model (Koster and Suarez, 1995), a mesoscale model, and 
a cloud ensemble model (Tao et al., 1996). This NASA Technical Memorandum documents 
the detailed treatment of the solar radiative transfer in the model. 

2. SPECTRAL BANDS 

Table 1 shows the spectral bands of the solar radiation model. There are eight bands in 
the ultraviolet and visible region (v>14280 cm 1 ) and three bands in the infrared region 
(v< 14280 cm 1 ). Also shown in the table are the absorbers and scatterers included in the 
calculation of solar fluxes in each band, as well as the equations used for computing the 
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Table 1. Gaseous absorption and Rayleigh scattering in the broadband radiation 
parameterizations and the equations used for transmittance and flux calculations. 
Clouds and aerosols ate included in all bands. 


Band 

Spectral Range 

Absorber/ 

Scatterer 

Equation 


(cm 1 ) 

(pm) 



1 

(44440-57140) 

(0.175-0.225) 

0 3 

Rayleigh 

3.12 

2 

(40820-44440 

35700-38460) 

(0.225-0.245 

0.260-0.280) 

0 3 

Rayleigh 

3.12 

3 

(38460-40820) 

(0.245-0.260) 

0 3 

Rayleigh 

3.12 

4 

(33900-35700) 

(0.280-0.295) 

0 3 

Rayleigh 

3.12 

5 

(32260-33900) 

(0.295-0.310) 

0 3 

Rayleigh 

3.12 

6 

(31250-32260) 

(0.310-0.320) 

0 3 

Rayleigh 

3.12 

7 

(25000-31250) 

(0.320-0.400) 

0 3 

Rayleigh 

3.12 

8 

(14280-25000) 

(0.400-0.700) 

0 3 

H 2 0 

Rayleigh 

3.12 

3.12 

9 

(8200-14280) 

(0.70-1.22) 

H 2 0 

0 3 * 

Rayleigh 

3.11 

3.12 

10 

(4400-8200) 

(1.22-2.27) 

H 2 0 

Rayleigh 

3.11 

11 

(1000-4400) 

(2.27-10.0) 

h 2 o 

3.11 


Total Spectrum 


0 2 

3.16 


Total Spectrum 


co 2 

3.19; 3.24 


* 0 3 absorption is folded into Band 8 
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transmission function. The first eight bands involves the 0 3 absorption and Rayleigh 
scattering. Band 8 also includes the weak absorption due to water vapor. Single values of the 
0 3 absorption coefficient, the water vapor absorption coefficient, and Rayleigh scattering 
coefficient are used in each of the 8 bands. Bands 9-11 include the water vapor absorption 
and Rayleigh scattering. Water vapor absorption in these bands is significant, and the k- 
distribution method is used. Band 9 also includes the weak 0 3 absorption. The 0 3 absorption 
in this band is folded into the absorption in Band 8. The absorption due to 0 2 and C0 2 is of 
secondary importance but occur in wide spectral ranges. Different approaches which 
compute only the reduction in fluxes are used. Clouds and aerosols are included in all bands. 

3. GASEOUS ABSORPTION 

3.1 High spectral-resolution calculations 

The absorption coefficients of water vapor, 0 2 , and C0 2 vary rapidly with wavenumber. 
Parameterizations for the absorption due to these gases are based on line-by-line calculations 
using the latest (1996) version of the molecular line parameters compiled at AFGL (Rothman, 
1987). The shape of an absorption line is assumed to follow the Voigt function. The Voigt 
half-width is a function of pressure, temperature and wavenumber. Because the shape of the 
far wings of a line is not well known, the absorption coefficient of a molecular line is set to 
zero at wavenumbers >10 cm 1 from the line center. The spectral resolution of the line-by- 
line calculations is 0.002 cm 1 , which is adequate for solar flux calculations. The 0 3 
absorption coefficient varies by several orders of magnitude in the UV and visible regions but 
rather smoothly with wavenumber. Absorption due to 0 3 is computed with a resolution of 
~0.003 pm. The absorption coefficient is interpolated from the spectral values given in WMO 
(1985). The effect of pressure and temperature on the 0 3 absorption is small and, hence, 
neglected. 

3.2 The k- distribution method 

The width of an absorption line of water vapor, 0 2 , or C0 2 is much smaller than the 
average spacing between absorption lines. It decreases linearly with decreasing pressure and 
is only ~0.05 cm 1 in the lower troposphere. As a result, the absorption coefficient varies 
rapidly within narrow spectral intervals, and accurate radiative transfer calculations require a 
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high spectral resolution. For the case of water vapor absorption, the high resolution line-by- 
line method requires flux calculations at >10 6 spectral points. For an atmospheric layer where 
temperature 0 and pressure p can be assumed constant, the wavenumbers with the same 
absorption coefficient are radiatively identical and can be treated as one entity. Fluxes at 
those wavenumbers with the same absorption coefficient need to be computed only once. 
Therefore, flux calculations can be greatly simplified by grouping the wavenumbers with the 
same absorption coefficient. Within a small interval 8v where the spectral variation of 
insolation (incoming solar radiation) is small, the integration of fluxes over wavenumbers can 
be replaced by that over the absorption coefficient. The mean transmission function of the 
interval 8v can be written as 

o° n 

T(w ' ) = J e~ kv(p ' e)H ' dv / <5v = J e~ k "'g(k; p, Q)dk - d)e~ ky (3.1) 

Sv 0 '=1 

where W is the absorber amount, k v is the absorption coefficient at the wavenumber V, g is 
the £ -distribution density function (cf. Arking and Grossman, 1972) such that the fraction of 
spectrum with the absorption coefficient between k — dk / 2 and k + dk / 2 is g(k)dk , and a } 
[=g(kj)Akj ] is the £ -distribution function. Chou and Lee (1996) has shown that the k - 
distribution method requires only 10 values of k (i.e. n= 10) for accurate calculations of 
solar heating due to water vapor using the £ -distribution method, instead of >10 6 spectral 
points using the line-by-line method. 

In the atmosphere where pressure and temperature change with height, the 
wavenumbers with a common absorption coefficient at a given height will not necessarily 
have a common absorption coefficient at other heights. These wavenumbers are no longer 
radiatively identical, and the k -distribution method cannot be used without applying 
assumptions. Pressure affects the width of a molecular line. Near the center of a line, the 
absorption is strong, and flux calculations are not sensitive to errors in the absorption 
coefficient. Distant from line centers, the absorption coefficient varies nearly linearly with 
pressure. On the other hand, temperature affects absorption primarily through its effect on 
the line intensity, which is constant for a given line. Thus, the temperature effect on 
absorption is rather smooth with respect to wavenumber. It has been shown by Chou (1986) 
that fluxes can be calculated accurately using the one-parameter scaling for the absorption 
coefficient, 
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k v (p,G) = k v (p r ,d r ) 


r P Y 


mo r ) 


(3.2) 


\PrJ 

where p r is the reference pressure, O r is the reference temperature, m< 1, and f(6,6 r ) is 
the temperature scaling function. As explained in Chou (1986), the optical path between the 
top of the atmosphere and the middle and lower troposphere is large near the center of 
absorption lines. It follows that solar fluxes in the middle and lower troposphere are 
primarily attributable to the spectral intervals between absorption lines, and the wing- 
approximation of (3.2) can be applied to accurately compute fluxes. In the stratosphere, the 
solar heating due to water vapor and O, is small, and underestimation of the absorption 
coefficient induced by the use of (3.2) does not have a significant effect on the heating rate. 
In the upper troposphere, the absorption is neither strong nor weak, and the error induced by 
the scaling is expected to have a large impact on flux calculations. However, if we choose p r 
to be the upper tropospheric pressure, then p / p r ~l in the upper troposphere, and the 
scaling introduces only a small error in the absorption coefficient. It has been found (Chou, 
1986) that the solar fluxes can be computed accurately by choosing p r = 300 hPa and 
m = 0.8. Flux calculations are not sensitivity to m for 0.7<m<0.9. The effect of temperature 
on the absorption of solar radiation is weak. For the absorption due to water vapor, it is 
approximated by 


f(6,d r ) = 1 + 0.00135(0-0,) (3.3) 

where 0,.=24O K. For the absorption due to 0 2 , and C0 2 , the temperature effect is neglected, 
i.e. f(0,O r ) = 1. 


Using the scaling of (3.2), Equation (3.1) becomes 

T(w) = ^a i (p r ,0 r )e~ kiH 

i = 1 


where 


w(p,0) = w' 


( X 

— rne r ) 

{Pr) 


(3.4) 


(3.5) 


Thus, the scaling of the absorption coefficient is reduced to the scaling of the absorber 
amount. The mean transmission given in (3.4) is independent of wavenumber, and flux 
calculations are greatly simplified. It is noted that the physical basis for the £ -distribution 
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method with the one-parameter pressure and temperature scaling is the same as for the 
correlated A' -distribution method (Wang and Ryan, 1983; Goody et al., 1989; Lacis and 
Oinas, 1991; Fu and Liou, 1992; Mlawer et al., 1997), i.e., the dominant effect of line wings 
on the absorption. The correlated A -distribution method correlates the values of A at 
different pressure and temperature through table look-up. The A -distribution method with 
the one-parameter scaling correlates values of A through (3.2). The accuracies of these two 
approaches are comparable, but the latter is much simpler than the former. 

3.3 Absorption due to water vapor in the infrared (Bands 9-11) 

Figure la shows the spectral distribution of the water vapor transmission function of the 
entire atmospheric column in the vertical direction averaged over 10 cm 4 , The atmosphere 
used is typical of mid-latitude summer conditions (Anderson et al., 1986). It shows that water 
vapor absorption varies strongly with wavenumber. The A -distribution method is used to 
parameterize water vapor absorption. For a spectral band with a width Av, the direct solar 
flux at a given pressure level p can be written as 

F(w)= je~ Kw S y dv (3.6) 

Av 

where S v is the extraterrestrial solar flux density at the wavenumber v multiplied by the 
cosine of the solar zenith angle, A is the absorption coefficient at p r and 0 r , and w is the 
scaled absorber amount, given by (3.5), above the pressure level p and in the direction of the 
solar beam. The flux-weighted mean transmission of a band is given by 

T(w)= \e~ Kw S y dv/ (3.7) 

Av Av 

If we divide a spectral band into m small intervals and apply the A -distribution method 
of (3.4), Equation (3.6) becomes 

m n m 

F(w) = ^Sj je~ k '"dv =X I>,,AA y ; e~ k ‘ W (3-8) 

j = 1 A xj i=l U=! 

where Sj is the mean flux density at the top of the atmosphere in the interval A V j, a , j is the 
A -distribution function in the spectral interval AVy, n is number of A -intervals, and 
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Wavenumber (cm ’) 


Wavenumber (cm ’) 


Figure 1. Spectral distribution of the transmission function of the entire column of a mid-latitude summer 
atmosphere in the vertical direction averaged over 10 cm 1 for water vapor, 0 3 , 0 2 , and CO,. 


m 

Av=^Av, (3.9) 

j = i 

By defining 

m 

h i = / S o (3-10) 

j = i 

which is the flux- weighted k-distribution function at p r = 300 hPa and d r =240 K, Equation 
(3.8) becomes 

F(w) = S 0 J d h f e~ kiW (3.11) 

r=i 

where S a is the total insolation at the top of the atmosphere (=extraterrestrial solar flux 
multiplied by the cosine of the solar zenith angle). Equation (3.11) is applied to compute the 
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water vapor transmission function in the infrared (Bands 9-11). Fluxes are computed for 
each kj, and the total flux is the sum of these fluxes weighted by the £ -distribution function. 
Table 2 gives the values of k and the k -distribution function used in the three IR bands. It is 
noted that if n is smaller than 10, say, 7, the calculated heating rate profile will oscillate 
strongly with height, as shown in Figure 1 of Chou and Lee (1996). It is also noted that the 
reason for dividing the infrared into three bands is the large change with wavenumber in the 
cloud and aerosol optical thickness and single-scattering albedo. 

Table 2. The flux-weighted k -distribution function of water vapor, h , in the infrared for p = 300 hPa and 
0 ; .=24O K. Sum of all It values is 0.52926. The unit of the absorption coefficient, k, is g 1 cm 2 . 


k -Interval 

k 

h (Band 9) 

h (Band 10) 

h (Band 11) 

1 

0.0010 

0.20673 

0.08236 

0.01074 

2 

0.0133 

0.03497 

0.01157 

0.00360 

3 

0.0422 

0.03011 

0.01133 

0.00411 

4 

0.1334 

0.02260 

0.01143 

0.00421 

5 

0.4217 

0.01336 

0.01240 

0.00389 

6 

1.3340 

0.00696 

0.01258 

0.00326 

7 

5.6230 

0.00441 

0.01381 

0.00499 

8 

31.620 

0.00115 

0.00650 

0.00465 

9 

177.80 

0.00026 

0.00244 

0.00245 

10 

1000.0 

0.00000 

0.00094 

0.00145 


3.4 Absorption due to water vapor in the PAR region (Band 8) 

The absorption due to water vapor in the photosynthetically active radiation (PAR) 
region between 14280 cm 1 and 25000 cm 1 (Band 8) is weak (Figure la), and the mean 
transmission defined by (3.7) can be approximated by 

T(w) * J (1 - k v w)S v dv / J S v d v = 1 - kw~ e~~ u (3.12) 

Av Av 

where k is an effective absorption coefficient, and w is the scaled water vapor amount given 
by (3.5). Using the line-by-line method, the mean transmission function T(w) given by 
(3.7) is first computed at 300 hPa and 240 K for a large range of w. For each w, an 
effective absorption coefficient is then computed by inverting (3.12), 

k(w) = — —In T(w) (3.13) 

w 
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The value of k is found to range from 0.00065 to 0.00080 g 1 cm 2 , and a mean value of 
0.00075 g 1 cm 2 is used to compute the water vapor absorption in Band 8 from (3.12). 

3.5 Absorption due to 0 3 in the UV and PAR regions (Bands 1-8) 

The eight bands in UV and PAR spectral regions as shown in Table 1 are relatively 
narrow. The range of the 0 3 absorption coefficient in each band is small and the mean 
transmission function can also be approximated by (3.12). By dividing the UV and PAR 
regions into 127 narrow bands with a width of ~0.003 pm and using the 0 3 absorption 
coefficient given in WMO (1986), the mean 0 3 transmission function is computed from (3.7). 
The same procedures described above for the determination of the effective water vapor 
absorption coefficient is applied to derive the effective 0 3 absorption coefficient, which is 
given in Table 3. It is noted that the effect of temperature and pressure on the 0 3 absorption 


Table 3. The spectral range, fractional solar flux, effective ozone and water vapor absorption coefficients, 
and Rayleigh scattering coefficient. Sum of S / S a is 1 . 



Band 

S / S a 0 3 Absorption 

H 2 0 Absorption 

Ray. Scattering 



Coefficient 

Coefficient 

Coefficient 



(cm-atrn) 1 

(g 1 cm 2 ) 

(hPa 1 ) 


uv-c 



1 

0.00057 

30.47 


0.00604 


2 

0.00367 

187.24 


0.00170 


3 

0.00083 

301.92 


0.00222 

UV-B 

4 

0.00417 

42.83 


0.00132 


5 

0.00600 

7.09 


0.00107 


6 

0.00556 

1.25 


0.00091 

UV-A 

7 

0.05913 

0.0345 


0.00055 

PAR 







8 

0.39081 

0.0572 

0.00075 

0.00012 

Near IR 

9 

0.32055 


See Table 2 

0.0000156 


10 

0.16536 


See Table 2 

0.0000018 


11 

0.04335 


See Table 2 

0.0 


* Enhanced by the 0 3 absorption in Band 9. 
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is small, and the one-parameter scaling is not applied to the 0 3 amount. Also shown in the 
table is the ratio of the extraterrestrial solar radiation of each band to that of the total 
spectrum, S / S a . Equation (3.12) is applied to compute the transmission due to 0 3 in Bands 
1 - 8 . 

3.6 Absorption due to 0 3 in the infrared 

As shown in Figure lb, the absorption due to ozone in the near infrared (Band 9) is 
restricted to a narrow spectral region next to Band 8 where absorption due to water vapor is 
weak. However, the absorption due to water vapor in the rest of Band 9 is not necessarily 
weak. Thus, calculations of the ozone heating in the near infrared would require adding 
another band in the broadband radiation model. Because the absorption due to water vapor 
and 0 3 is weak in Band 8, the 0 3 absorption in Band 9 can be folded into the absorption in 
Band 8 as if it were the absorption due to another absorber. To do so, the effective 0 3 
absorption coefficient of Band 8, k% , is replaced by 

k = k $ + Ak (3.14) 

where A k satisfies 

(l - e~ Akn ' ) J S v dv = J(l -e~ kvW )Sydv (3.15) 

Av 8 Av 9 

so that the absorption of solar radiation due to 0 3 in the near infrared is correctly computed, 
where W is the ozone amount, and Av g and Av 9 are the widths of Bands 8 and 9, 
respectively. For a wide range of W found in the atmosphere, the value of A k ranges 
between 0.0032 and 0.0033 (cm-atm) STP \ Therefore, the ozone absorption coefficient in 
Band 8 is enhanced by 0.0033 (cm-atm) STP 1 to take into account the absorption by ozone in 
the infrared. The 0 3 absorption coefficient of Band 8 given in Table 3 includes this extra 
absorption. This approach to computing the ozone absorption in the near infrared requires 
no extra computing time. 

3.7 Absorption due to 0 2 

The absorption due to 0 2 occurs in narrow spectral intervals, but is not necessary weak 
near band centers (Figure lc). Using the line-by-line method, Chou (1990) computed the 
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mean 0 2 transmission function in the A and B bands centered at 13150 cm 1 and 14510 cm 
and fit the transmission function by a square-root approximation, 


T(w) = e~ l4w 


(3.16) 


The absorption in the weak bands centered at 7890 cm 1 and 15870 cm 1 was not 
included. To include the absorption in all those bands, the mean transmission function of 0 2 
at p r = 300 hPa and 0,.=24O K in the spectral regions 7600-8050, 12850-13190, 14310- 
14590, and 15730-15930 cm 1 is computed from (3.7) using the line-by-line method, and the 
effective mean absorption coefficient k is computed from 


k(w ) = - 



ln7Xw) 


(3.17) 


The value of k varies between 0.000135 (cm-atm) STP 1/2 and 0.000155 (cm-atm) STP 1/2 for a 
large range of the scaled 0 2 amount, w, encountered in the atmosphere. A mean value of 
0.000145 (cm-atm) SXP 1/2 is used in the broadband radiation model. 

The 0 2 absorption bands are located between water vapor absorption bands. Thus, 
water vapor has little effect on the solar heating due to 0 2 . Instead of computing fluxes in the 
0 2 bands, the reduction in fluxes due to 0 2 absorption is computed from 

AF(w) = s(l-e - ^) (3.18) 

The insolation, S, in the spectral regions 7600-8050, 12850-13190, 14310-14590, and 
15730-15930 cm 1 is 6.33% of the total solar flux at the top of the atmosphere. 

The 0 2 transmission function, Equation (3.16), is similar to that used by Kiehl and 
Yamanouchi (1985), except the pressure scaling of (3.5) is different. It is band-averaged, 
does not follow Beer’s law, and, hence, cannot be applied to the 8 -Eddington approximation 
for multiple-scattering calculations. We apply (3.16) only to clear atmospheres. For cloudy 
atmospheres, a simple scaling approach is used to compute the 0 2 heating, which is addressed 
in Section 6.4. 
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3.8 Absorption due to C0 2 

The absorption due to C0 2 is small but overlaps significantly with the absorption due to 
water vapor in certain spectral regions (Figures la and Id). By assuming that C0 2 heating in 
the stratosphere is not important, Chou (1990) computed the reduction of fluxes due to C0 2 
based on the pressure scaling of (3.5) and ignoring C0 2 absorption in the weak bands at 6400 
cm 1 and 7000 cm 1 and in the strong 2300 cm 1 band, where the insolation is weak. A 
deficiency of this approach is that the heating rate in the upper stratosphere was greatly 
underestimated due to the use of the pressure scaling. To demonstrate the importance of the 
C0 2 heating in the stratosphere, Figure 2 shows the line-by-line calculated heating rate profile 
(solid curve) in the spectral region 2200-2400 cm 1 and 3500-3760 cm' 1 , where C0 2 
absorption is strong. The heating rate increases from 1 °C day 1 at 0.1 hPa to ~ 4.5 °C day 1 at 
0.01 hPa . For an atmospheric model which extends to the middle atmosphere, it is desirable 
to include this large heating. 



Figure 2. The effect of C0 2 absorption in the spectral 
regions 2200-2400 cm 1 and 3500-3760 cm 1 on the 
solar heating. Overlapping with water vapor 
absorption is included in the calculations of C0 2 effect. 
The solid curve is line-by-line calculations, and the 
dashed curve is the parameterization using the 
precomputed table of AF(p c / /J Q ,p ) given by 
(3.19). The C0 2 concentration is set to 350 pprnv, aid 
the solar zenith angle is set to 60°. 


In the 2200-2400 cm 1 band, water vapor absorption is negligible; whereas in the 3500- 
3760 cm 1 band, water vapor absorption is very strong in the troposphere but weak in the 
stratosphere, and the effect of C0 2 on the tropospheric heating can be neglected. It follows 
that the absorption due to C0 2 in these spectral regions is not sensitive to water vapor amount, 
and an average water vapor path can be used in the parameterization for the flux reduction 
due to C0 2 . By using a fixed water vapor path, the C0 2 -induced flux reduction at a given 
pressure level becomes a function only of the C0 2 concentration and the solar zenith angle, 
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(3.19) 


A F(p c / p a ,p ) = J e~^ ip) [\ - e~ rApJp °’ p) ]s v dv 


where 


1 p 

r w (p) = =-\k v (p)q(p)dp 
pg J 0 

(3.20) 

*c(Pc/ Fo’P)= Pc \YvWP 

II a J 

(3.21) 


k and y are the absorption coefficients of water vapor and C0 2 , p c is the C0 2 concentration 
independent of space, q is the water vapor mixing ratio of a typical midlatitude summer 
atmosphere, fl Q is the cosine of the solar zenith angle, and p is taken to be 0.5. Using the 
line-by-line method, a table for A F(p c / p o ,p) normalized by S f) is precomputed for 
p c / p o ranging from 250 pprnv to 3300 ppmv and p ranging from 0.01 hPa to 1000 hPa. 
The C0 2 heating in 2200-2400 cm 1 and 3500-3760 cm 1 derived from the precomputed table 
of A F(p c / p o ,p) is shown as the dashed curve of Figure 2. It agrees very well with the line- 
line-line calculations (solid curve). 

For computing the C0 2 heating outside the strong bands at 2200-2400 cm 1 and 3500- 
3760 cm 1 , the approach of Chou (1990) is followed. Similar to (3.6) and (3.11), fluxes in 
the C0 2 bands are given by 

n m 

F(w,u ) = J e - {k ' w+y ' u) S r dv = S 0 XX /t './~ a ' H+7 '" ) (3.22) 

'=1 ;= i 

where w and u are the scaled water vapor and C0 2 amounts, and /?, y is the flux-weighted k- 
distribution function with the absorption coefficients kj and y j for the two absorbers, 
respectively. Equation (3.22) requires ( mxn ) sets of flux calculations, which becomes 
computationally very expensive when used in weather and climate models. Other approaches 
than the k -distribution method have to be used for efficient flux calculations when 
overlapping of absorption is involved. 

The reduction in transmission due to C0 2 absorption can be expressed as 

AT(w,u) = J e~ kvH (l - e~ 7rH )S v dv / J S v dv (3.23) 
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Since AT(w,u ) is a function only of w and u, it is pre-computed using the line-by-line 
method for large ranges of w and u encountered in the atmosphere. The integration is over 
the total solar spectrum. The absorption due to C0 2 in 2200-2400 cm 1 and 3500-3760 cm 1 , 
is excluded in computing AT(w,u). Fluxes in Band 10 are then simply computed from 

F(w,u ) = F(w)- S a AT(w,u), (3.24) 

where F(w) is computed from (3.11). 

4. MOLECULAR AND PARTICULATE SCATTERING 
4.1 Rayleigh scattering 

Rayleigh scattering decreases with decreasing wavenumber according to ~v 4 . Spectral 
values of Rayleigh scattering extinction coefficient are taken from a World Meteorological 
Organization report on atmospheric ozone (WMO, 1986). For each of Bands 1-10, a mean 
effective Rayleigh scattering coefficient is computed according to the same procedures as in 
the computation of the mean effective water vapor absorption coefficient in the PAR region, 
(3.13). The Rayleigh scattering coefficients for Bands 1-10 are given in Table 3. 

4.2. Cloud and rain single-scattering properties 

Radiation routines for atmospheric models usually use only a few broad spectral bands, 
with effective optical properties parameterized for each band. The single-scattering co- 
albedo, 1 — Q) v , varies by several orders of magnitude in the solar spectrum. It is, therefore, 
difficult to derive an effective (0 for a broad spectral band that can be applied to both strong 
and weak absorption situations. There are a number of approaches to deriving (O . Slingo 
and Schrecker (1982) derived (0 by weighting (O v with the solar insolation and linearly 
averaging over a spectral band. Espinoza and Flarshvardhan (1996) used a square-root 
approximation to derive (0. Fu (1996) used a mix of linear and logarithmic averaging 
depending on the strength of absorption in a given spectral band. The absorption of solar 
radiation by clouds depends not only on (O v but also on water vapor and cloud ice/water 
amounts. It is clear that no optimal method has been found for deriving (0 over a broad 
band. It can only be empirically determined based on the amounts of cloud particles and 
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water vapor encountered in the Earth's atmosphere. Let us define the linear and logarithmic 


averaging over a band as 


(1 - (O ' ) = Xd - 0) v )P v S v Av / ^P V S V Av 

(4.1) 

log(l - a>" ) = X 1°S(1 - co v )P v S v A v / X £ A Av 

(4.2) 


where S v is the solar insolation at the top of the atmosphere (TOA), and Av is a narrow 
spectral interval where optical properties can be treated as constants. The effective mean 
single-scattering co-albedo of a band is then computed from 

(1- (O) = 77(1 -*»') + (1 “ f?)(l-G>") (4.3) 

where T) is the weight ranging from 0 to 1. The weight Tj should be close to 1 for the 
spectral bands with weak absorption and should decrease as absorption increases. Chou et al. 
(1997) determined empirically the weight in the various spectral bands that minimize errors 
in flux calculations. 

Compared to Q) v , the extinction coefficient ( ) and the asymmetry factor ( g v ) vary 
rather smoothly with V. Their effective mean values over a wide spectral band can be 
accurately approximated by 


0 = 2/5AAv/£SvAv 

(4.4) 

g = / £ 01 Ji, S v Av 

(4.5) 


Theoretical considerations and radiative transfer calculations have shown that cloud 
single-scattering properties are not significantly affected by details of the particle size 
distribution and can be adequately parameterized as functions of the effective particle size 
(Fu, 1996; Hu and Stamnes, 1993; Tsay et al., 1989). Following Slingo and Schrecker 
(1982), we parameterize the single-scattering properties for a broad spectral band by 


P = a 0 +ajr e 

(4.6) 

1 ~(0 = b o +b 1 r e + h 2 r e 2 

(4.7) 
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(4.8) 


g = c a + c 1 r e + c 2 r e 


2 


where a , b , and c are regression coefficients, and r e is the effective particle size defined to 
be proportional to the ratio of the total volume of cloud particles to the total cross-sectional 
area, A c of cloud particles. For spherical water droplets, the effective size is given by 


= J r 3 n(r)dr / J r 2 n(r)dr 


3 C 
4 P M . 4 


(4.9) 


where p w is the density of water. For hexagonal ice crystals randomly oriented in space, the 
effective size is shown by Fu (1996) to be 


• _ 2 V3 C 

r ' " 3p, 4 


(4.10) 


where p, is the density of ice, and C is the cloud water/ice mass concentration per unit 
volume. The cloud optical thickness, T c , is then given by 

T c =j3Cz (4.11) 

where z is the geometric thickness of a cloud layer. 

The size, shape, and refractive indices are different for ice crystals and water droplets. 
The extinction coefficient of ice clouds is smaller than that of water clouds because ice 
crystals are much larger than water droplets. The single-scattering albedo and asymmetry 
factor of ice clouds and water clouds are also different. Therefore, we parameterize the 
single-scattering properties separately for ice and water clouds. 

The spectral data (O v , (5 V , and g v calculated by Fu (1996) for ice clouds and by Tsay 
et al. (1989) for water clouds are used to derive ft), (5 , and g from (4.1)-(4.5) as functions of 
phase and size. By assuming hexagonal ice crystals randomly oriented in space, Fu (1996) 
computed the single-scattering parameters of ice clouds using the improved ray-tracing 
method of Yang and Liou (1995). A total of 28 size distributions derived from in situ 
aircraft measurements were used in the calculations, which included samples from the FIRE 
(First ISCCP Regional Experiment) and CEPEX (Central Equatorial Pacific Experiment) field 
campaigns. The mean effective size of ice crystals, r ( , ranges from 20 pm to 130 pm. For 
water clouds, the single-scattering parameters of 5 water clouds were computed by Tsay et al. 
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(1989) from the Mie theory assuming spherical droplets and a log-normal size distribution. 
The effective radius of water droplets, r w , ranges from 4 pm to 20 pm. 

Figures 3-5 show the results of regression using (4.6)-(4.8). The coefficients a , b, and 
c are shown in Tables 4-6, respectively. It can be seen in Figure 3 that the extinction 
coefficient varies weakly with spectral band but strongly with particle size. For large particles, 
it is independent of wavelength for both water and ice clouds. The results shown in Figure 3 
for ice clouds are indistinguishable among the four bands. Due to a large particle size, the 
extinction coefficient of ice clouds is significantly smaller than that of water clouds. It is 
interesting to note that in spite of a large difference in the particle size distribution, the ice 
particle extinction coefficients for the 28 clouds nearly fall onto a single curve of (4.6). 
Figure 4 shows the single- scattering co-albedo of Bands 9-11. The single-scattering co- 

-5 

albedo of Band 1-8 is very small (<10 ) and is set to 0. For a given particle size, the single- 
scattering co-albedo varies by 3 orders of magnitude among the three IR bands. For a given 
spectral band, it is ~ 8 times larger for ice clouds than for water clouds. The asymmetry 
factor, shown in Figure 5, varies between 0.78 and 0.94 for different bands and particle sizes. 
Since the shape of size distribution has little effect on the single-scattering properties (Flu and 
Stamnes, 1993; Slingo and Schrecker, 1982), the results shown in Figures 3-5 should be 
representative for a wide range of clouds with various particle size distributions. 

For rain drops, a single drop-size distribution is used to compute the spectral 
distribution of the single-scattering properties. Fu et al. (1995) used a truncated constant- 
slope gamma function (Manton and Cotton, 1977) with a minimum size of 60 pm, a 
maximum size of 1800 pm, and a mean effective size of 162 pm to compute the single- 
scattering properties. The mean effective extinction coefficient, single-scattering co-albedo, 
and asymmetry factor for the 1 1 spectral bands are shown in Tables 4-6. They are applied to 
all other drop-size distributions. The results are provided to us by Prof. Qiang Fu 

4.3. Aerosol single-scattering properties 

The aerosol optical thickness, single-scattering albedo, and asymmetry factor are not 
parameterized. They need to be specified as functions of height and spectral band for input 
to the radiation routine. 
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Table 4. Coefficients a o and of tlie parameterization, Equation (4.6), for the cloud extinction 
coefficient. The units of a o and <7, are, respectively, m 2 g 1 and m 2 g 1 pm. 


Spectral Band 

Ice Cloud 

Water Cloud 

Rain 


<*o 

<h 

a Q 

0\ 

a o 

<h 

1-8 

3.33e-4 

2.52 

-6.59e-3 

1.65 

3.07 e-3 

0.0 

9 

3.33e-4 

2.52 

-1. Ole-2 

1.72 

3.07 e-3 

0.0 

10 

3.33e-4 

2.52 

-1 ,66e-2 

1.85 

3.07 e-3 

0.0 

11 

3.33e-4 

2.52 

-3.39e-2 

2.16 

3.07 e-3 

0.0 



Figure 3. The extinction coefficient of ice and water clouds as a function of the effective particle size, 
r, and r w are, respectively, the effective size of ice and water particles. The data points aie 
derived from Equation (4.4), whereas the curves are the regression using Equation (4.6). 
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Table 5. Coefficients b 0 , b h and b 2 of the parameterization, Equation (4.7), for the cloud single-scattering co-albedo. 
The units of bi and b 2 are, respectively, pm' 1 and pm' 2 . 
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7.49e-l 1.20e-3 -3.67e-6 7.94e-l 8.32e-3 -2.33e-4 0.891 

7.61e-l 1.42e-3 -3.96e-6 7.45e-l 1.37e-2 -3.82e-4 0.948 

8.41e-l 1.26e-3 -3.85e-6 8.35e-l 2.57e-3 5.52e-5 0.971 




Figure 4. The single-scattering co-albedo of ice and water clouds as a function of the effective particle 
size, r, and r w are, respectively, the effective size of ice and water particles. The data points aie 
derived from Equations. (4. l)-(4.3), whereas the curves are the regression using Equation (4.7). 
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Figure 5. The asymmetry factor of ice and water clouds as a function of the effective particle size. 7) 
and r are, respectively, the effective size of ice and water particles. The data points are derived 
from Equation (4.5), whereas the curves are the regression using Equation (4.8). 
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5. SURFACE ALBEDO 


The surface albedo is an input parameter and is separately specified for the UV and 
PAR region (Bands 1-8) and the infrared (Bands 9-11). It is also separately specified for 
direct and diffuse fluxes. Therefore, a set of four surface albedos must be specified as input 
to the radiation routine. 

6. FLUX CALCULATIONS 

As mentioned in the previous sections, single values of the gaseous absorption 
coefficients and particle single-scattering properties are used for each of the 11 bands, except 
for the water vapor absorption coefficient in the IR region (Bands 9-11) where the k- 
distribution method is used. Therefore, one set of fluxes in the atmosphere and at the surface 
are computed for each band in the UV and visible regions (Bands 1-8), and 10 sets of fluxes 
corresponding to 10 values of the water vapor absorption coefficient are computed for each 
of the IR bands (Bands 9-11). Reflectivities and transmissivities of each atmospheric layer 
are first computed using the 8-Eddington approximation (Joseph et al., 1976), fluxes are then 
computed by using a two-streant adding method for a composite of layers. The total flux at 
each pressure level and at the surface is the weighted sum of these fluxes, 

38 

F(p) = ^Vi F i(P) ( 6 - 1 ) 

i=i 

where y/ j is the fraction of extraterrestrial solar flux contained either in one of the first 8 
bands (S / S 0 given in Table 3) or in one of the 30 k-intervals in the three IR bands ( h given 
in Table 2). The sum of I//) is 1. 

6.1 Layer optical properties 

For each atmospheric layer and spectral band, the effective optical thickness, single- 
scattering albedo, and asymmetry factor are computed from 


I*. 

i 

(6.2) 

/X T , 

i i 

(6.3) 


(6.4) 


i i 
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where the summation is over all gases and particles, / . They include ozone, water vapor, 
clouds, aerosols, and atmospheric gases (Rayleigh scattering). The single-scattering albedo is 
set to 1 for cloud scattering in the UV and visible regions (Bands 1-8) and for Rayleigh 
scattering in all bands, whereas it is set to 0 for ozone and water vapor. Cloud ice and water 
particles are allowed to mix in a layer. Absorption due to C0 2 and 0 2 is treated differently 
and is addressed in Section 6.4 

6.2 Layer reflectivities and transmissivities 

The transmissivity T and reflectivity R of a layer illuminated by a direct beam are 
computed from the S-Eddington approximation as a function of p D and the optical properties 
T, 0) , and g defined by (6.2)-(6.4). King and Harshvardhan (1986) have shown that errors 
in transmissivity and reflectivity induced by this approximation are ~5%. For diffuse 
insolation, the transmissivity T and reflectivity R are also computed from the 8-Eddington 
approximation but incident angles are approximated by a single value of 53°, 


R ~ R(/J ) 

(6.5) 

T~T(p ) 

(6.6) 


where JJ=c os(53°). 

To evaluate the errors induced by this approximation, we use the discrete-ordinate 
multiple-scattering algorithm of Stamnes et al. (1988) to compute the transmissivity and 
reflectivity of a layer for various incident angles. The diffuse transmissivity and reflectivity 
are derived by averaging the direct-beam transmissivity and reflectivity over 271 solid angles, 


l 


2 J R{p)pdp 
0 

(6.7) 

1 

2 J T(p. pd/J 

(6.8) 


o 


where jJ is the cosine of the incident angle depart from the vertical. Errors induced by (6.5) 
and (6.6) are shown in Figure 6 for large ranges of the optical thickness and single-scattering 
co-albedo. Figure 6 shows that the approximation introduces an error < 0.02 in both the 
reflectance and transmittance for wide ranges of the optical thickness and the single-scattering 


22 



Log 10 ( 1 


albedo, and an asymmetry factor of 0.8. For other values of the asymmetry factor, results are 
similar. 


Reflectance Difference Transmittance Difference 




Figure 6. The left panel is the diffuse reflectivity computed from (6.5) minus that computed from (6.7). 
The right panel is the diffuse transmissivity computed from (6.6) minus that computed from 
(6.8). The cloud asymmetry factor is set to 0.85. 


6.3 Two-stream adding approximation and flux calculations 

Let us consider a layer divided into an upper sublayer a and a lower sublayer b . The 
direct-beam transmissivity T , the direct-beam reflectivity R , the diffuse transmissivity T , and 
the diffuse reflectivity R of each of the two sublayers are computed from the 5-Eddington 
approximation. By separating the direct and diffuse components of the radiation, the 
reflectivity and transmissivity of the layer when illuminated by direct radiation can be 
computed from (Chou, 1992) 

p? / (6.9) 

(1 -R a R b ) 

^ ) = ^ ^ (ju 0 ) + 7^ fc {e- ) + [^ (^^ ) - ]} / (6.10) 

(i ~R a R b ) 
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where T{jA 0 ), e r/ ' Ll " , and T(fl a ) — e z/ “° are the total, direct, and diffuse transmissivities, 
respectively. The subscripts a and b denote the upper layer and the lower layer. 

The reflectivity of the layer when illuminated by diffuse radiation from above, R ah , is 
generally different from that when illuminated from below, R ab . They can be computed 
from 


R ab =R c ,+ T a R b T a /(i-R a R b ) 


R ab- R b +T b R a T fc/C 1 R a R b ) 


( 6 . 11 ) 

( 6 . 12 ) 


where the overbar and double-overbar denote illumination from above and below the layer, 
respectively. 

T 

When illuminated by a direct beam radiation at the top of the layer, the reflectivity, R , 
and the transmissivity, T ^ , at the interface of the sublayers a and b are given by 

R r ={ e - T ^R b ^ 0 ) + [T a (ji 0 )-e- T ^p b }l{\-R a R b ) ( 6 . 13 ) 

T l =e- Tj ^ +{e- Tj ^R a R b ^ 0 ) + [T c ^ 0 )-e- Tj ^ ]}/(\-R a R b ) ( 6 . 14 ) 

It follows from (6.13) and (6.14) that the upward and downward fluxes at the level 
/ + 1 / 2 (the lower boundary of the layer i in the atmosphere) can be computed from 


-h.z/t'o 


R 


i+l,k 


/ 


(6.15) 


W T = s{e‘ v "" K/ +u M ) + [r,, (#<„) — «■' 

_ R 1 ,i R /+U-) 

F Mn l = Se-" jl >‘- +s{e-"‘ l “-R Ul R iH ^ 0 ) + [T lJ ^ a )-e-"‘ l ^l (6.16) 

0- — R l,i R i+l,k) 


where S is the insolation at the top of the atmosphere, the index k denotes the surface, and 
the subscripts (1, i) and ( / + 1 , k) denote the regions above and below the level / + 1/2, 
respectively. The transmission 7) ( and reflection R lj of the region above the level are 
computed by adding layers from the top of the atmosphere to the layer i using (6.10) and 
(6.12) according to the procedures given in Lacis and Hansen (1974). Similarly, the 
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reflectivities R j+l k and R )+1 k are computed by adding layers from the surface up to the 
layer i + 1 using (6.9) and (6.11). 


Fluxes are computed for each of the 8 bands in the UV and PAR regions and for each 
of the 30 £ -intervals in the IR bands. Total fluxes are the sum of these fluxes given by (6.1). 


6.4 Flux reduction due to C0 2 and 0 2 

The reduction of the downward flux due to 0 2 and C0 2 in clear sky conditions, AF'Kir, 
is computed from (3.18), (3.19), and (3.24). To simplify flux calculations in cloudy 
conditions, the effect of 0 2 and C0 2 on solar heating is approximated by 


A fKm 
AFKm 


(p) = AF l c , r (p), 


(p) = AFKirip) 


F" e, dd(p) 

F" e 'c,r ( P ) 


for p< p c 
for p > p c 


(6.17) 

(6.18) 


where F net is the total net downward flux without 0 2 and C0 2 absorption, and p c is the cloud- 
top pressure. 

7. CLOUD OVERLAPPING 

Clouds could occur at various heights with different fractional cover and optical 
properties. Nearly all radiative transfer algorithms used in atmospheric models apply only to 
a plane-parallel (i.e. horizontally homogeneous) atmosphere. Horizontal inhomogeneity is 
not allowed. A straightforward approach to dealing with a partial cloudiness situation is to 
divide a region into sections. Within each section, an atmospheric layer is either free of 
clouds or filled totally with a homogeneous cloud. Radiative fluxes are then computed for 
each section, and the total solar heating is the sum of all sections weighted by their fractional 
cover. Depending upon the number of cloud layers and the way these clouds overlap, 
computational costs of this approach could be huge. To simplify the computations, a cloud 
that partially fills a layer is usually smeared over the entire layer, and the optical thickness, T, 
is adjusted by a factor dependent on the fractional cloud cover. As demonstrated by a 
number of studies (e.g. Harshvardhan and Randall 1985; Cahalan et al., 1994), the effect of 
T on radiation is highly nonlinear, and it is not appropriate to scale T (or cloud water path) 
linearly by cloud cover. In the NCAR (National Center for Atmospheric Research) 
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community climate model CCM, the optical thickness is scaled by f 1 ’ 5 , where / is the 
fractional cloud cover (Briegleb, 1992). In the Goddard GCM, Sud et al. (1993) scaled the 
optical thickness by inverting a simple reflection function for clouds when illuminated by 
diffuse radiation and assuming no absorption by clouds. As can be expected, the effective 
optical thickness corresponding to a cloud smeared to cover the entire layer is a function of 
T, /, and the solar zenith angle, as well as whether scaling is based on reflection or 
absorption. Scaling of T is further complicated by the ways that clouds in various layers 
overlap. 

7.1. Maximum-random overlapping 

Treatments of cloud overlapping in climate studies differ from model to model. In our 
radiation parameterization, a combination of maximum and random cloud overlapping 
scheme similar to that of Geleyn and Hollingsworth (1979) is used (Chou et al., 1997). 
Clouds are identified as high, middle, or low separated roughly by the 400-hPa and 700-hPa 
levels. Cloud layers close to each other are likely to be related; therefore, cloud layers in each 
of the three height groups are assumed to be maximally overlapped (left panel of Figure 7), 
whereas overlapping among the three groups is assumed to be random. Scaling of X is 
applied only to the maximally-overlapped clouds within each of the three height groups, but 
not to the randomly-overlapped clouds among different height groups (right panel of Figure 
7). Within each height group, clouds are smeared over the extent of the maximum cloud 
amount, f m , of that height group by scaling T (right panel of Figure 7). By assuming 
random overlapping of clouds among different height groups, the atmosphere is then divided 
into 2" sections, where n < 3 is the number of height groups containing clouds (right panel 
of Figure 7). Within each section, a layer is either totally cloud-filled or cloud-free. Fluxes 
are first computed for each section and then summed over all sections weighted by the 
fractional cover of the individual sections. 

As examples, the fractional cover of the clear section, A c/r , is computed from 


A 


v clr 




(7.1) 
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where the subscripts 1, 2, and 3 denotes the upper, middle, and lower cloud groups, 
respectively. Similarly, the fractional cover of the section with only high clouds is computed 
from 




(7.2) 


Cloud Overlapping Optical-thickness Scaling 



Figure 7. The maximum-random cloud-overlapping scheme. Clouds are identified as high, middle, and low 
separated roughly by the 400 liPa and 700 hPa levels. Clouds are assumed maximally overlapped 
within each height group and randomly overlapped among different height groups (left panel). 
Scaling of the optical thickness applies only to the maximally-overlapped cloud layers (right panel). 
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7.2. Scaling of optical thickness 

Cloud optical-thickness scaling used to “smear” the clouds to f m within each height 
group is based on the maximum-random cloud overlapping assumption shown in Figure 7. 
To simplify the scaling of T, we assume that reflection is more important to climate than 
absorption, since the former affects the total heating of the earth-atmosphere system, but the 
latter only redistributes the heating within the atmosphere, as well as between the atmosphere 
and the surface. As stated in Sections 6.2 and 6.3, our radiation algorithm requires 
calculations of reflection/transmission of a cloud layer separately for direct and diffuse 
radiation. By smearing cloud cover over the extent of the maximum cloud cover of 
respective height-group, the optical thickness of a cloud layer is scaled by factors % s and % s . 


II 

u 

for direct radiation 

(7.3) 


for diffuse radiation 

(7.4) 

albedos of the layer remain the 

same before and after scaling, 


R(r b ,r e ,Ho) = fR(.T,r e ,iio) 

for direct radiation 

(7.5) 

R(t f ,r e ) = f'R(t,r e ) 

for diffuse radiation 

(7.6) 


where p o is the consine of the solar zenith angle, R is the albedo averaged over the entire 
solar spectrum for direct incident radiation, R is the albedo for diffuse incident radiation 
given by 

t 

R(T,r e ) = 2jR(T,r e ,ii)/ud/J (7.7) 

o 

f is the nomialized cloud cover given by 

f = flf m (7-8) 

and fJ is the cosine of the zenith angle of an incident beam. Using the scaling of (7.3) and 
(7.4), a layer with a fractional cloud cover / and an optical thickness T is reduced to a layer 
with a fractional cloud cover f m and an equivalent optical thickness T h for direct radiation 
and tj for diffuse radiation. It is noted that the scaling is a function of the incident angle, 
and that optical thickness is scaled separately for direct radiation and diffuse radiation. 
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Using values of r e representative of the water and ice cloud samples used in Section 4.2, 
CO and g are computed from (6.3) and (6.4) for the individual bands given in Table 1 . 
Cloud albedos R and R averaged over the solar spectrum are computed using the discrete- 
ordinate multiple- scattering algorithm of Stamnes et al. (1988). The scaling functions % s 
and X s are then derived from (7.3)-(7.6) by giving /' , R( T,r e ,fA o )and R( t,r e ). Sample 
calculations show that X s and X s depend only weakly on r e and the phase of cloud 
particles. Therefore, they are reduced to X.Af i UAU) and ^(/’,T) by taking the mean 
values of the scaling functions for various particle sizes and phases. 

The above scaling ensures that the reflection of a cloud layer is the same before and 
after the scaling. However, errors are expected in the calculations of atmospheric and surface 
heating when there are other cloud layers present. In addition to the solar zenith angle, cloud 
cover, and optical thickness, an optimal cloud scaling should, in principle, also depend on the 
fractional cover and optical thickness of other cloud layers. This situation makes scaling of 
the optical thickness extremely difficult. To simplify the problem, we made simple 
adjustments to the scaling functions, Xs and Xs> so that errors in the fluxes at the TOA and 
the surface are minimized for the simplest cloud situation, which has two cloud layers in only 
one of the three height-groups. For the two-cloud situation, scaling is applied to the cloud 
with a smaller cloud cover (see right panel of Figure 7). It is noted that scaling is not 
required when there is only one cloud layer in a given height-group. The scaling functions 
are adjusted according to 


X = Xs + (f-Xs)S 

(7.9) 

X = Xs + (f-Xs)S’ 

(7.10) 


where 5 and d are constants empirically (trial and error) chosen to minimize errors in the 
fluxes at the TOA and at the surface. The errors are defined as the differences between the 
flux computed with the scaling of T and the flux computed with the sky divided into 
horizontally homogeneous sections (three sections in the two-cloud situation). The adjusted 
scaling functions preserve the shapes of Xs a °d Xs- The results for the scaling functions are 
saved as tables for later use. It is found that optimal values of 8 and 8 are 0.2 and 0.3, 
respectively. The results for the scaling functions X and X are shown in Figure 8. It can be 
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seen that X(f , T,p o ) < f , and X f as T — >0. For an optically thick cloud, X and X 
are not sensitive to T . 


The scaling functions derived for the two-cloud situation is applied universally to all 
situations even if there are more than two cloud layers in a given height group. It implicitly 
assumes that, in a given height group, accurate scaling of T is not critical if there are two 
cloud layers above. 

With the sky divided into 2" sections, fluxes are computed for each section. Instead of 
(6.1), the total flux is computed from 


38 


F(P) = ^¥i 


i=i 


;'=i 


A F U 


(p) 


(7.11) 




Figure 8. Scaling functions for the cloud optical-thickness, X an d X, f° r direct (shown only for p Q = 
0.3) and diffuse radiation, respectively. The value of f is the normalized cloud cover defined 

j., i 5 

by (7.8). Numbers in the figure are the cloud optical thickness. The dots are X = J use d m 
the NCAR CCM model (Briegleb, 1992). 

where N = 2” , and A ■ is the fractional cover of tire j th section. Tire sums of A- s and 
iff, s are both 1. It needs to compute the reflectivities and transmissivities of each of tire 
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atmospheric layers only once but not N times. Depending upon the number of sections N , 
the total computing time is about three tunes larger than the computing time for a single 
section. 

8. SAMPLE CALCULATIONS 

8.1. Effect of absorption and scattering in minor bands 

Concerned with the huge time required for computing radiative transfer in weather and 
climate models, solar heating in minor absorption bands has often been neglected. The 
individual contributions of these minor bands to the atmospheric heating is small, but 
collectively they are not negligible (Chou, 1999). The minor absorption and scattering which 
are often neglected include the absorption due to water vapor in the photo synthetically active 
radiation (PAR) spectral region (Band 8) and in the thermal infrared (1000-2600 cm 1 ), the 
ozone absorption and Rayleigh scattering in the near infrared (Band 9), and 0 2 and C0 2 
absorption in wide spectral ranges. Ramaswamy and Freidenreich (1991) have shown the 
importance of the two minor water vapor regions. Table 7 shows the effects of these minor 
absorption and scattering bands on the surface solar heating for the typical midlatitude 
summer atmosphere with a C0 2 concentration of 350 ppmv, a solar zenith angle of 60°, and a 
surface albedo of 0. Using a spectral resolution of 0.002 cm 1 for water vapor, 0 2 , and C0 2 , 
and 0.003 pm for 0 3 , the total effect of these minor bands is computed with the discrete- 
ordinary scattering algorithm. The surface solar heating is reduced by 17 W m . The 
atmospheric solar heating increases by a similar amount, ~15 W m (not shown) when these 
minor bands are included. This amounts to about 10% of the total atmospheric heating, 
which is not negligible in the solar heating calculations. The results of the parameterizations 
for these minor bands are also shown in the table, which are accurate compared to the detailed 
calculations. 

8.2. Clear-sky fluxes and heating rates 

Clear-sky fluxes and heating rates computed with the parameterization are validated 
with high spectral-resolution calculations. The results are shown in Table 8 for fluxes and 
Figure 9 for the heating rate. Calculations are for the mid-latitude summer atmosphere with a 
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Table 7. Effects of minor absorption and scattering on the downward surface solar flux from detailed 

high spectral-resolution calculations and parameterizations. Fluxes are computed for a typical 
mid-latitude summer atmosphere. The C0 2 concentration is set to 350 ppmv. The solar 
zenith angle and the surface albedo are set to 60° and zero, respectively. Units are W m 2 . 



Detailed Calculations 

Parameterization 

Water vapor in PAR 

-2.53 

-2.42 

Ozone in near infrared" 

-0.56 

-0.54 

Oxygen in near infrared" 

-4.29 

-4.18 

C0 2 in middle infrared" 

-3.30 

-3.57 

Rayleigh in near infrared 

-3.12 

-2.82 

H 2 0 in thermal infrared " 

-3.09 

* * * 

Total 

-16.89 



Overlapping with water vapor absorption is included. 

1000 cm L <v<2600 cm 1 . 

The effect cannot be assessed because a separate parameterization excluding 
the thermal infrared is not done. 


C0 2 concentration of 350 pprav. The solar zenith angle is set to 60°, and the surface albedo is 
set to 0.2. For the high resolution calculations, the spectral interval is chosen to be 0.003 pm 
for 0 3 and 0.002 cm 1 for water vapor, 0 2 , and C0 2 . Only the absorption due to water vapor, 
0 3 , 0 2 , and C0 2 are included. Scattering due to clouds, aerosols, and gases are excluded. 

Absorption of solar radiation in the atmosphere (ATM) and net downward fluxes at the 
top of the atmosphere (TOA) and the surface (SFC) are shown in Table 8 for different bands 
and the total spectrum. The parameterization introduces an error of <1% in the TOA and 
surface fluxes and the atmospheric heating. Figure 9 shows the heating rate as functions of 
the logarithm of pressure (left panel) and of pressure (right panel). The parameterization 
computes very accurately the heating rate below the 1-hPa level but with an error of 5% of 
the peak heating above. Most of the error is attributable to the use of the pressure scaling in 
the calculations of absorption due to water vapor, 0 2 , and C0 2 . 
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Table 8. Net downward solar fluxes at the top of the atmosphere (TOA) and at the surface (SFC) and the 
radiation absorbed hi the atmosphere (ATM) computed using high spectral resolution calculations 
and the parameterization. The atmosphere used is typical of a mid-latitude summer with a C0 2 
concentration of 350 ppmv. The solar zenith angle and the surface albedo are set to 60° and 0.2, 
respectively. Calculations include the absorption due to water vapor, 0 3 , C0 2 and 0 2 covering the 
spectral range from 0.175 pm to 10 pm. Scattering is excluded. Units are W m 2 . 



High Spectral-Resolution 

Parameterization 


TOA 

SFC 

ATM 

TOA 

SFC 

ATM 

Band 1-7 

47.7 

33.7 

14.0 

46.7 

33.0 

13.7 

Band 8 

217.7 

203.2 

14.5 

218.0 

203.9 

14.1 

Band 9 

186.7 

137.8 

48.9 

186.5 

134.9 

51.6 

Band 10 

101.5 

51.0 

50.5 

100.3 

51.7 

48.6 

Band 11 

28.0 

7.8 

20.2 

27.6 

8.2 

19.4 

Total 

581.5 

433.3 

148.2 

579.1 

431.7 

147.4 



Heating Rate (C/day) Heating Rate (C/day) 


Figure 9. The heating rate profiles of a mid-latitude summer atmosphere computed with the line-by- 
line method (solid curves ) and the parameterization (dashed curves). The heating is expressed as 
a function of the logarithm of pressure in the left panel and as a function of pressure in the right 
panel. The solar zenith angle is set to 60°, surface albedo is set to 0.2, and the C0 2 
concentration is set to 350 ppmv. Absorption due to water vapor, 0 3 , 0 2 , and C0 2 ate included. 
Scattering is excluded. 


33 





8.3. Cloudy-sky fluxes and heating rates 


Fluxes and heating rates in a cloudy atmosphere derived from the parameterization and 
the detailed absorption and scattering calculations are compared. A deck of stratus cloud 
covering the entire sky is inserted between 800 hPa and 920 hPa. The effective radius of 
liquid cloud particles, r w , is set to 12 pm, and the optical thickness in the 0.55 pm visible 
region, T vjs , is set to 9.7. Based on results from Mie calculations for those r w and T vjs , the 
extinction coefficient in the visible region, /3 rjs , is 0.1303 g m 2 , and the total cloud water 
amount, Cz given by (4.11), is 74.4 g m 2 . The cloud layer is divided into 5 sublayers with a 
thickness of 24 hPa. For each sublayer, the cloud water amount is 14.9 g m 2 , and the optical 
thickness in the visible region is 1.94. 


For the detailed calculations, the line-by-line method and the discrete-ordinate 
scattering algorithm are used. The spectral distributions of the cloud extinction coefficient, 
/3 V , the single-scattering albedo, (O v , and the asymmetry factor, g v , for r w = 12 pm are 
derived from Mie calculations. The spectral distribution of the optical thickness is scaled 
according to 


r = t 

V VIS 


A 

vA is J 


(8.1) 


For the parameterization, the band-averaged cloud extinction coefficient, the single-scattering 
albedo, and the asymmetry factor are derived from (4.6)-(4.8), and the optical thickness is 
derived from (4.11). Fluxes are computed using the two-stream adding approximation given 
in Section 6.3. 

Table 9 shows the net downward fluxes at the TOA and surface and the solar heating of 
the mid-latitude summer atmosphere using the detailed calculations and the parameterization. 
The solar zenith angle and the surface albedo are set to 60° and 0.2, respectively. The 
discrepancy between the two sets of calculations is very small in the surface flux and is ~2% 
in the total TOA flux. Figure 10 shows the heating rate profiles of the two sets of 
calculations. They agree very well both inside and outside the cloud. 
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Pressure (hPa) 


Table 9. Same as Table 8 except scattering by gases and cloud is included. A deck of stratus cloud 
with an optical thickness of 9.7 in the visible region is located between 800 and 920 hPa. 
Units are W in . 


Detailed Calculations Parameterization 



TOA 

SFC 

ATM 

TOA 

SFC 

ATM 

Band 1-8 

138.7 

105.8 

32.9 

137.8 

107.3 

30.5 

Band 9 

110.6 

62.9 

47.7 

106.4 

61.9 

44.6 

Band 10 

78.5 

17.9 

60.6 

74.1 

18.5 

55.7 

Band 11 

27.1 

0.7 

26.4 

28.0 

0.1 

27.9 

Total 

354.9 

187.2 

167.7 

346.4 

186.9 

159.6 



Figure 10. Same as Figure 9, except a 
cloud deck is inserted between 800 and 
920 hPa. The cloud optical thickness is 
9.7 in the visible region, and the 
effective radius of cloud liquid chops is 
12 pm. 
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